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. _ , Abstract. A two-dimensional circular quantum billiard with unusual boundary 

conditions introduced by Berry and Dennis (J Phys A 41 (2008) 135203, Ref. [2]) 
is considered in detail. It is demonstrated that most of its eigenfunctions are strongly 
localized and the corresponding eigenvalues are close to eigenvalues of the circular 
billiard with Neumann boundary conditions. Deviations from strong localization are 
also discussed. These results agree well with numerical calculations. 
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1. Introduction 

> , 

G^ , It is common to consider a quantum problem as integrable (but not necessarily 
i separable) or chaotic, depending solely on the properties of its classical counterpart. For 
^ I example, the quantum circular billiard is integrable because classical mechanics inside 
^ I a circle is integrable. In such an approach one does not even specify the exact form 
0\ ', of boundary conditions, which are indispensable for the existence of discrete spectrum 
^ ' in quantum billiards. Often this approach is correct; for example a circular billiard 
remains integrable for both Dirichlet and Neumann boundary conditions. Nevertheless, 
^1 it is known that this is not always the case. A notable example is a rectangular billiard 
5^ I with Dirichlet conditions imposed on part of the boundary, and Neumann conditions on 
the complimentary part (e.g. and references therein). Though the rectangular billiard 
is classically integrable, this version of the quantum problem is neither integrable nor 
chaotic and, in fact, has many features in common with classically pseudointegrable 
systems. In particular, in [1] it was shown that the eigenfunctions of such systems 
have a strong resemblance with integrable eigenfunctions whose quantization gives the 
positions of the energy levels with reasonably good precision. 

Recently, a different system with similar properties was introduced in [2]. The 
problem consists in finding eigenfunctions and eigenvalues for the Helmholtz equation 
in polar coordinates [2] 

(A + A;2)^(r,0) = (1.1) 

with Robin boundary conditions (mixed boundary conditions) on the circle of radius R: 

|-vl/(r,0)|,=^ = A/(0)v^(i?,0) . (1.2) 
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When the boundary function /(0) is a constant, the problem remains integrable 
(e.g. [3]). New phenomena appear when 

• f{(f)) is a smooth function of the polar angle (f), and 

• the prefactor A is not a constant but is proportional to the momentum k. 
In [2] the case 



and A = k was briefly considered. Fig. 7 of that paper shows an eigenfunction with 
energy E = 97.206986712 which looks like the (slightly shifted) standard Bessel function 
with azimuthal quantum number 8 and radial quantum number 6. The 6**^ zero of Jg{x) 
is xo = 27.8892694, corresponding to energy Eq = x^/S = 97.226418653. Compared 
with the unit mean distance between levels, the difference between these is rather small: 
E-Eo = -0.01943. 

The purpose of this paper is to show that this not a coincidence. We argue that 
almost all energy levels of the problem (11.11) . (II. 2p . and (I1.3P are close to eigenvalues of 
the circular billiard with Neumann boundary conditions 



The corresponding wavefunctions look, roughly speaking, like the usual Bessel function 
solutions of (II. ip but slightly shifted with respect to the circle centre. 

The main reason for such behaviour is the strong localization of wave functions. 
Namely, the recurrence relation determining the eigenfunction coefficients (see equation 
(12. 2p below), is equivalent to the one-dimensional discrete Schrodinger equation with a 
pseudo-random distribution of on-site energies. If these energies were truly random, this 
problem would correspond to Anderson localization in one dimension O [6l [7] for which 
all states are localized. Though the energy distribution is deterministic in our problem, 
it has strong pseudo-random properties, and many of the predictions of localization 
theory remain valid. 

The plan of the paper is the following. In Section [2] the formal solution of the 
problem in terms of recurrent relations derived in [2| is discussed. Properties of 
the pseudo-random energy distribution and its relation with the Lloyd model of the 
product of random matrices are investigated, and the localization length is discussed. 
Properties of strongly localized states and the construction of a local perturbation theory 
are treated in Section [3l States which are locahzed close to the boundaries of the 
allowed region differ from strongly localized states. In Section HI states with large 
azimuthal number are shown to be adequately described by a continuous semiclassical 
approximation, whose approximate quantization condition is derived in Section [51 As for 
any dynamical model, the construction of the semiclassical trace formula is of interest, 
and this is done for our model in Section El We conclude in Section [71 Details of longer 
calculations are given in the Appendices. 



/(0) = COS0 



(1.3) 



J'^ix) = . 



(1.4) 
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2. Localization 

The solution of the above problem may be written as the following formal series [2] 



where x = kR and cosm0 (resp. sinm0) are chosen for states symmetric 
(resp. antisymmetric) with respect to the symmetry transformation —0. 

The boundary conditions fll.2p and fll.3l) are fulfilled provided the coefficients 
obey the recurrence relations [2j 

(2.2) 

where 

for all m = 0, 1, 2, . . .. The initial values may be chosen as follows [2] 

ao = 1, ai = po(a;) , (2.4) 

for symmetric functions and 

ao = 0, ai = 1 (2.5) 

for antisymmetric ones. 

For integrable models such as the circular billiard with Neumann boundary 
conditions, only one term in the series (12.11) is non-zero. For chaotic problems like the 
stadium billiard, all coefficients are non-zero, and in the mean the numbers am/Jm{x) 
can be considered as independent Gaussian random variables [Ij. Below we demonstrate 
that for our problem, almost all states are localized. This means that one coefficient am 
is much larger than all the others which, roughly speaking, decrease exponentially from 
the centre of localization. 

It is evident that the recursion relations (12. 2p define the discrete Schrodinger 
equation, and can be rewritten in the form of the transfer matrix 




(2.6) 



In the semiclassical approximation x — > cxd and m < x, pm can be approximated by 
using the standard asymptotics of the Bessel functions |8] 



^{x^-m^r'/\os^Ux) , (2.7) 



where 



Thus 



^mix) = V x"^ — w? — m arccos . (2-8) 

X 4 



Pn^x) ~ -\/ 1 - tan($^(a;)) . (2.9) 
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If, for all m, ^m{x) are independent random variables distributed uniformly between 
and TT, the variables Pm{x) are independent random variables with the Cauchy 
distribution 

w 

P(P) = ( 2^ 2^ ' (2-10) 

with 



w = ^l-m'^/x'^ . (2.11) 

In such a case, Eq. (12.61) determines the soluble Lloyd model [9] of the product of random 
matrices, for which it is known that for almost all initial conditions the Lyapunov 
exponent of the product, defined 

A = hm ^< + , (2.12) 

m— >oo Tfl 



is non-zero, and for the distribution (I2.10p its value is given by [TO] 

A = In (w + Vl + w2) . (2.13) 

In parallel with the Lloyd model we also consider a closely related deterministic model 
(sometimes called the Maryland model [H]) defined by the same recursion relations 
(12. 2p . but with the function 

= sin(^) tan(m^) (2.14) 

with a certain constant 9. When 6' is a 'good' irrational multiple of vr, this model should 
have properties close to the Lloyd model. In particular, the Lyapunov exponent is given 
by (12.131) . These two models are compared in Fig. [1^). On the other hand, when 
is rational (or irrational but with a good rational approximation), the situation differs 
considerably from a random model (see e.g. |12j). 

The Lyapunov exponent (I2.12p is a self-averaged quantity only in the limit m — oo. 
For a finite sample length (i.e. large but finite m < L), it is a random variable with a 
certain distribution. In [11] it was shown that this distribution in the Lloyd model is 
Gaussian with mean A and variance approximately given by 

- ^ . (2.15) 

In Fig. [TJo the distribution of the Lyapunov exponents for the Lloyd model with L = 100 
is shown. For a given realization of a random sequence, the Lyapunov exponent with 
finite L has a certain value, but the mean value and variance over many realizations 
agree well with fl213D and fl2J[5D . 

When boundary conditions for the Schrodinger equation are imposed at large 
distances, almost all eigenstates, are localized and for one-dimensional systems the 
localization length / equals the inverse of the Lyapunov exponent [15] 



l = j. (2.16) 
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Figure 1. (a) Logarithm of modulus of tlie coefficients versus m for the 
deterministic model ()2.14p with 6 — 7r\/5/4 (the upper solid line) and for one realization 
of the Lloyd model (|2.10p with w — 2 sin 6 (the lower solid line) . For clarity the upper 
curve is shifted by 5 units. The dashed line is the Lloyd model prediction y — Am where 
A is the Lyapunov exponent given by (|2.13p . (b) Distribution of the Lyapunov exponent 
for 100 realizations of the Lloyd model with the length L = 100 and parameters the 
same as in (a). The dashed line is the Gaussian distribution with parameters (|2.13p 
and (jTTS)) 



For the Cauchy distribution with the width defined by (12. lip , the local localization 
length is 



1 



l{t) = ln(v^r^ + V2^) (2.17) 



with t = m/x where m is the center of localization and x is the eigen-momentum. For 
a large range of m/x, the localization length is close to 1 (e.g. for m = .5x, I = 1.28) 
but near m = a; it diverges. 

Of course, this statement is valid only for the pure random Lloyd model with the 
width as in (12. lip . In our case, the phases ^m{x) in (12.81) are not random, but rather 
quickly varying functions of m. First of all, the oscillatory asymptotics of the Bessel 
functions (12. 7p are valid only when 

|m| < X . (2.18) 

Therefore, wave functions with eigen-momentum x (if any exist) can be localized only 
in the interval 

< m < a; . (2.19) 

Secondly, within this interval, the best one can expect is that the phases (12. 8p are 

pseudo-random, provided that their derivative over m is not a rational multiple of tt, 

nor too close to one. This means that in regions where 

d , . m M , , 

— $^ X = - arccos — = — tt 2.20 
dm X N 

one cannot expect good localization of eigenfunctions. 
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Figure 2. Natural logarithm of coefficient modulus In \am\ versus m for a symmetric 
eigenvalue with x « 32.50302694 (black circles) and for x ~ 32.50524735 corresponding 
to a zero of Jg{x) (open squares). Dashed lines have the slopes given by expression 
([2l7l) with m = 9. 



Nevertheless, for most values of m, the phases ( 12.8P are pseudo-random when 
considered modulo vr (cf. [I3]), and it is natural to assume that wave functions obeying 
(12. 2p are localized (at least in a certain interval of m). 

In Fig. O the absolute values of the coefficients am corresponding to an eigenvalue 
are presented on the logarithmic scale. It is clear that they correspond to an eigenstate 
localized at m = 9 and decaying exponentially from this point with the localization 
length close to the one given by Eq. (12.171) . When the momentum is even slightly 
different from a true eigenvalue, the coefficients grow exponentially from this point. 



3. Strongly localized states 

We will call states localized at a certain point m*, not too close to the boundary and 
X, strongly localized states. As the localization length / is of order 1, in general, these 
states consist of one large component am* and all other components have to be small: 

am ~ e-l"*l/' . (3.1) 

Therefore it is natural to develop a perturbation series on the number of amplitudes in 
the vicinity of m*. 

We assume that all coefficients except the one with m = m* are zero, and am* = 1- 
Then from the recurrence relations (12. 2p . it follows that the eigenvalue x* has to be a 
zero of the derivative of the Bessel function with m = m*, 

JUx*) = 0; (3.2) 
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The next approximation consists of taking into account terms with m = m* ± 1, which 
leads to the 3x3 equation 

/ -2p„._i(a;) 1 \ \ 

1 -2pm*ix) 1 1 =0. (3.3) 

\ 1 -2p^.+i(x) J ya+i J 

Of course, this and the following equations can be solved numerically, but it is more 
convenient to perform some calculations beforehand. We shall see that solutions of any 
of these equations lead to x close to x* from (13. 2p . Therefore, one has to know the values 



of Jm*+k(x*) and JL 



m*+k 



'x*). From standard recursion relations for the Bessel functions 



(3.4) 



zJ'^{z) ±nJn{z) = ±zJn^i{z), 
it follows that the Bessel functions can be written in the form 

Jn+k{z) = Jn{z)Rk,n{z) - J n-l{z) Rk^l ,n+l{z) (3.5) 

where Rk,n{z) represents a certain polynomial of degree k in 1/z called Lommel's 
polynomial [8]. In particular, when x is a zero of J'^{x), direct calculations give 

X X 

m{m + 1) 



and 



x^ 



XI 



<-^m('^) 



-1 + 



m(m — 1) 



x^ 



(3.6) 



Jm+2(yX) 



m-2 



X] 



2m(m + 1) 



x^ 



2(m + 1) 



- IJ , Jrn-2{X) = Jm{x) 

2m{m + l)(m + 2) 



2m{m — 1) 



x^ 



X 



X'^ 



(3.7) 



2(m ^ 1) _|_ 2m(m — l)(m — 2) 



X 



Using the properties of Lommel's polynomials or directly applying asymptotic formulas 
for the Bessel functions (as in Appendix C ), one finds that the first terms of semi classical 
expansion of Jp+k, calculated at a zero x of J^, are the following 

1 



Jp+k{x) = Jp{x) Tk{u) + 



2x(l-u^] 



[kuTkiu) + ik\l - u^) - l)f/fc_i(M)] +0(x~2)(3.8) 



and 



^ [kuUk-i{u) + eTk{u)]^+0{:. 



,-2\ 



(3.9) 



where u = p/x and Tk{u) and Uk{u) are the Chebyshev polynomials of the first and the 
second kind respectively, 

sin {k + 1)9 



Tkicos 



cos k9 , ?7i,(cos 



sin 6 



(3.10) 
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From these formulas it follows that 

Pm*+kix*) ^-Z^ sin 6* tan /c6' (3-11) 

with cos^ = m*/x*, and the last expression is an odd function of k: pm*+k{x*) = 
—Pm*~k{x*). In Appendix A it is shown this property implies that all (2p + 1) x (2p + 1) 
determinants as in (13. 3p vanish in the semiclassical limit, implying that their zeros are 
always close to the zero of the central element. 

For determinants of small size, the first order correction 

x = x* + — (3.12) 

X* 

can be calculated analytically. For 3x3 determinant in (13. 3p one obtains that 

= "(l-w2)(2«4-4«2 + 3) ^^-^^^ 
and for the 5x5 determinant 

4(8^^^ - 12u^ + 5m' + 1) 

^~ (1 - m2) (64^8 - 288^6 + 500^4 - 388^^ +115) ^ ' ^ 

where 

m* 

u = — . (3.15) 

X 

For example, for the eigenvalue x = 32.50302694 represented in Fig. [2], one has the 
following chain of approximations: 

• For the determinant 1x1, i.e. for the zero of ^9(^1) = one has 

X - xi ^ -0.0022 . (3.16) 

• For the solutions of 3 x 3, 5 x 5, and 7x7 determinants one finds numerically the 
following approximations 

X-X3 ^ -0.0013 , z-xs ^ -0.00022 , Z-X7 ^ -0.000038 .(3.17) 

Such perturbation expansions cannot converge uniformly for two different reasons. 
First, although it is reasonable to expect that eigenvalues obtained from such small 
determinants give a good approximation to the true eigenvalue (cf. (I3.17P ) for strongly 
localized states, the parameter of such an expansion is a pure number of the order of 
e~^ where A is the inverse of the localization length. Due to considerable fluctuations of 
the latter, it is difficult to give a precise a priori bound of the accuracy of such series. 
Second, starting from the center of the localization, the correct boundary conditions - as 
in (12. 4p and (12. 5p - are not taken into consideration. Nevertheless, boundary corrections 
have to be of the order of e"^''^ where 6N is the distance of the center of the localization 
to the boundary and can often be ignored for strongly localized states. In particular, 
this leads to the almost degeneracy of even and odd states, which differ only by the 
boundary values (12.40 and (12.50 . which is well-confirmed by numerics. For example, the 
odd state corresponding to even state in Fig. [2] has the momentum x = 32.50302689 
and the difference between even and odd states equals 5.34 10"^. 



Near integrable systems 



9 




Figure 3. Graphical definition of notation for the Bessel function addition theorem 



The above perturbation approach also gives information about the corresponding 
wave functions. Although for strongly localized states, the coefficients am decay 
exponentially from the localization center, numerically this decrease is not so quick, 
and the influence of the first corrections are noticeable. 

We now estimate these corrections from the 3x3 matrix of (13. 3p . Simple 
calculations reveal that, in the semiclassical limit, the wave function has the form (in 
this approximation) 

^(r, 0) ~ J„,(A;r)e^"^'<^ + a {Jm'+i{kr)e'^""+'^^ - J^*_i(A;r)e'("^*-^)<^) , (3.18) 

where m* is the position of localization center and the correction a is given by 
_ 1 

"~ 2(l-(m*/x)2) ■ ^^-^^^ 
The fact that the coefficient in front of the Jm*-i term is the opposite of that of the 



Jm*+i term is not related to the approximation used. From [Appendix A| it follows that 
in the semiclassical limit, when (13. lip is fulfilled, this will always be the case. 

To see better the main effect of this correction, we consider the familiar Bessel 
function addition theorem [8]. According to this theorem, the Bessel function shifted 
by a vector e can be expanded as follows: 

oo 

JmHe'"^'^ = Yl ^n(e)J™+n(r)e^("^+")^ (3.20) 

n=— oo 

where (r, (p) and {w, ip) are the polar coordinates of a point with respect to the axis 
defined by the direction of the shift e (see Fig. [3]). 

Taking into account only the smallest of e terms, and using J_„(x) = (— l)"J„(x), 
one concludes that 

JUkw)e'^^ ^ Jo{ke) [J^(fcr)e''"<^ + P{Jm+i{kr)e'^'^^'^^ - J^-i(fcr)e^("^-^)<^)] (3.21) 
where 

/3 = ^ . (3.22) 






IMM*» S'lHIl 





Figure 4. Left: Wave function of the even state witli x sa 32.50303 sliown in 
Fig. [21 Riglit: Wave function '^{r, (p) = J<j{kr) cos9(/) with /ci? w 32.50525 satisfying to 
Neumann boundary conditions. 



Comparing this equation with (13.181) we conclude that the most noticeable effect of 
corrections is a small shift to the left of the state with m = m*, the value of this shift 
being of the order oil/k. For illustration we present in Fig. H] the numerically-computed 
wave function corresponding to the eigenstate shown in Fig. [2l 

4. Continuous approximation 

As the localization length grows when the localization centre m* is close to x, the 
structure of the wave functions in this region differs considerably from the that discussed 
in the preceding Sections. An example is shown in Fig. [51 The absence of a sharp 
localized peak together with the existence of a flat part between two vertical lines are a 
characteristic feature of such a state, clearly differing from the strongly localized state 
of Fig. [21 Nevertheless, its eigen- momentum x = 28.40353635 is still close to the lowest 
zero of ^26(^0), at xq = 28.4181, the difference being x — xq = —0.0145. 

To explain such behaviour, we use a semiclassical-type approximation based on 
the fact that, in the region m ^ x, the effective 'potential' in the discrete Schrodinger 
equation (Q, V{m) = 2p„( smooth function of m (see Fig. [6]). When m> x and 

X ^ 00, the Bessel function Jm{x) can be approximated by an asymptotic expression 



[8J similar to (12^1) : 




1 



(4.1) 



where 




(4.2) 
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Figure 5. CoefBcicnts for an even state with x — 28.40353635 (black circles connected 
by solid line). Vertical solid line indicates the position of eigen momentum x. The 
abscissa of the vertical dashed line equals y/2x. 



from which it follows that 



Pm{x) + ITT— 2 ^ • 4.3 

This formula agrees well with the exact Pm{x) for large m, except in a region close to 
m = X, where the following uniform approximation is useful (see e.g. [8J): 

where Ai{y) is the Airy function of argument 

y = —{m — x) (4.5) 

and the width v = (x/2)^/^. 

In particular, V{m) diverges for m = rUp, determined by the closest zero of the 
Bessel function, Jmp{x) = 0. For the choice of x in Fig. [6], nip ~ 22.79270. In the 
uniform approximation (14. 4p . 

mp = x- 2-^l^rix^'^ + O^x'^'^) (4.6) 

where rj is the modulus of the first zero of the Airy function, r] 2.338. 
We thus have the discrete Schrodinger equation 

V{m)am = a-m+i + a-m-i (4.7) 

with smooth 'potential' V{m) = 2pm{x). It is therefore natural to look for its particular 
solution am in a semiclassical form (see e.g. [T7] ) 

am = A(m)e^*('") (4.8) 
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Figure 6. 'Potential' 2pm{x), with pm{x) given by (12. 3|) calculated at a; = 28.41807 
as a function of m. The solid line is the semiclassical approximation (|4.3p . defined 
only when m > x. The dashed line indicates the value of energy in the continuous 
Schrodinger equation ()4.15|) . 



with certain smooth functions $(m) and A{m). 

From the considerations below it follows that $(m) — * 
and am-i may be approximated as the Taylor series 

a^+i ^ (A(m) + A'(m))e'l*('")+*'("^)+^*"('")l , 

and 

am-i ~ {A{m) - A'(m))e'[*("^)'*'('")+5'f"M] . 
Equating coefficients in the difference equation (H? 
$'(m) = arccospm(2;) 

and 

1 



oo when x ^ oo, and am+i 

(4.9) 



in this approximation 



A{m) 



(4.10) 
(4.11) 
(4.12) 



(1 - pi{^)y/' ■ 

The oscillating solutions with real $(m) exist when 

\pUx)\ < 1 (4.13) 

and values of m, where Pm{x) = ±1, play the role of turning points. Outside the interval 
fl4.13p one solution is growing and the other is decaying. 

Close to the right turning point where pm{x) ~ 1, the continuous Schrodinger 
equation may be constructed directly from the discrete one (14.71) by expanding am±i 
into a series over the derivatives with respect to m: 



(4.14) 



Near integrable systems 13 
One then gets the Schrodinger equation 

d2 



dm? 



+ 2-V{m)]ip{m) = (4.15) 



with 'wave function' ^/^(m) = am{x). 

The solutions (14. 8 p are exact analogues of the usual semiclassical solutions of this 
Schrodinger equation, and the matching formulas connecting the decaying solution with 
the oscillating one may be obtained in the same way as in the standard semiclassical 
case (see e.g. [T8]). 

We denote 

{arccospm(x) when |pm(x)| < 1 

In {\pUx)\ + VpU^) - l) when |p^(x)| > 1 ' ^^'^^^ 
and let 1712 be the right turning point i.e. the solution of the equation 

Pm,{x) = l. (4.17) 
In the semiclassical limit, m2 can be well-estimated from the approximation (14.30 

m2 = V2x-^. (4.18) 

When m > 1712, the semiclassical solution which tends to zero when m ^ 00 has the 
form 

^ exp (- I q{t)dt\ , (4.19) 



2{pl{x) - 1)V4 
and when m < m2, it can be approximated as 



1 / /"""^ tt' 

cos / q{t)dt - - ) . (4.20) 



When m > X and m is not too close to x, one can use semiclassical approximation (14.31) 
for pm{x)- For m > \/2x in the leading approximation, one obtains 

1 



2(n2 - 2)1/4 
and for x < m < \/2x 

1 



exp 



(2-^2)1/4 



COS 



-X j \n (yW^ + Vt^^) dt 



i arccos ( V — 1 ] dt — 

J u 



X 



4 



(4.21) 



(4.22) 



where u = m/x. As expected, the momentum x plays the role of 1/h. 
The left turning point, mi, is defined by the largest solution of 

p„,(x) = -l. (4.23) 

As this point is close to the pole (14. 6p . it follows that as a; ^ 00, 

mi ^nip + l (4.24) 
where is given by (14.61) . 
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Figure 7. (a) Coefficients for the state of Fig. [SI normalized to obey the condition 
(|4.19p . The sohd Hnes are semiclassical formulas (|4.19p and ()4.20p with x = xq. The 
vertical lines are the same as for Fig. \6\ (b) Wave function of the state shown in Fig. [S] 



Fig. shows the result of numerical integration of the semiclassical formulas 
f l4.19p and f l4.20p with the exact potential V^(m). The discontinuity close to ^/2x is 
due to the usual inapplicability of semiclassical formulas close to turning points. This 
figure also includes the values of coefficients am shown in Fig. [5] on the logarithmic scale. 
To compare them with semiclassical formulas, we normalize the coefficients so that 
their asymptotics agrees with (14.191) . To achieve this we multiply all a^ by a factor 
such that equals the value predicted by (14.191) . The agreement is good , confirming 
the applicability of the continuous approximation to describe states localized close to 
m = X. 

Though for the states we describe, many coefficients are non-zero, the 
corresponding wave functions are quite simple (cf. Fig. [7b))- The point is that almost 
all of these coefficients correspond to the Bessel functions Jm{x) with m > x, which 
decay exponentially inside the circle. 

The behaviour of for large m described in this Section is not, of course, specific 
to states localized with m close to x. Rather, it is evidently generic, and all states have 
a similar form for very large m. For example, the state in Fig. [2] shows clear exponential 
localization, but only in the finite interval m < x. In Fig. [8^), the same state is 
plotted over a larger interval. Starting from m = x, it deviates from pure exponential 
localization. From Fig. [HJo) it follows that its large-m behaviour is well described by the 
continuous approximation discussed above. 
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Figure 8. (a) Coefficients on the logarithmic scale for the same state with 
X « 32.50303 as in Fig. ^ but over a larger interval (black circles). Vertical solid 
and dashed lines indicate respectively m — x and m — x\f2. (b) The same as in (a) 
but in the usual scale and only for high m. 



5. Approximate quantization condition 

If our problem were a problem of a particle in a potential well, the semiclassical 
quantization condition determining the energy levels inside the well would take the 
form 

(5.1) 



/ g(m) dm b\ = ix-p 

J mi 4 



where Si is the phase shift associated with the left turning point and p is an integer. 

For the discrete equation (14.71) . there are two main differences from the standard 
case. First, there are no true local bound states. The exact position of the energy levels 
depends on the precise behaviour of Pm{x) far from the left turning point. Secondly, 
in the usual continuous Schrodinger equation, the integral of the momentum over the 
forbidden zone is real, but in our case, between and mi it also has an imaginary part 
equal to vr. 

Taking into account this additional phase, in this case the approximate quantization 
condition has the form 



rm2 

/ q{m) dm + n{mp} — 5 = up 

J mi 



(5.2) 



where nip is the position of the pole (14. 6p . {/} denotes the fractional part of /, and S 
denotes the sum of all phases, which we assume to be a slowly-varying function of the 
momentum. 

As the fractional part of a number differs from that number only by an integer, we 
can rewrite the above expression as 



m, 



rm2 

H — / g(m) dm — 5 = P 

71" J mi 



(5.3) 
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with integer P. 

The main term of this approximate quantization condition, when x — > cxd, takes 
the form 



x{l + -I) + g{x) -S = P, 

TT 

where / is the elhptic integral flB.161) 

1=1 arccos Vf^ — Idt . 



(5.4) 



(5.5) 



g{x) is a function which increases more slowly than x. A few of its low order terms can 
be read from flB.191) in Appendix B, 



— ^^/^ + — lnx + C(l) . 



24/3- '6,— — (5-6) 
The states localized in the region close to x are analogues of extreme whispering gallery 
states. For billiards with, say, Neumann boundary conditions, such states are quantized 
as X ~ P for integer P. The above discussion demonstrates that for the problem under 
consideration these states are strongly perturbed, and in the strong semiclassical limit 
they have momentum 

x^xP (5-7) 
with integer P. Here the factor 



X 



;i + _/)-! -0.9 . 



TT 



(5.^ 



A physical picture of these states follows. Consider the value of the wave function along 
the boundary of the circle. According to (12.11) , 

vl/(i?,0) = 5^a^e''-^ (5.9) 



Estimating coefficients from (14.221) . one finds 
^{R,(f)) ~ y^cos 



arccos 



X 



TT 



lldt-- 



(5.10) 



Calculating the sum in the saddle point approximation, we obtain that the saddle point, 
rusp, obeys the equation 



COS( 



m. 



sp 



1 . 



x^ 



Therefore the saddle point exists if 
cos > . 

When this condition is satisfied, the saddle point is the following 



sp 



x^ 



cos^ 



and a simple calculation shows that when — 7r/2 < </> < 7r/2, 



exp 



ix / \/l + cos^ fdf 
_ Jo 



(5.11) 
(5.12) 
(5.13) 

(5.14) 
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This wave function corresponds to a local evanescent mode propagating along the 
boundary 



exp 



ix [ ^/TT 
. Jo 



cos^ fdf — k{R — r) cos(0) 



(5.15) 



In general, such waves exist when the boundary function Af{(j)) entering the boundary 
conditions (11. 2p is positive [2]. 

Taking into account that at points ±vr/2, when cos0 = 0, the angular momentum 
of wave fl5.14p equals x, from continuity one concludes that when cos is negative the 
wave function along the boundary should have the form 

^(i?,0) ~ expix0 . (5.16) 

As the wave function has to be univalued, the quantization condition of such a state is 

1 + cos2 0d0 1 = 2nP (5.17) 



X TT 



7r/2 



with a certain integer P. From (1B.16P and (IB.lSp . it follows that this condition coincides 
with CT . 



6. Trace formula 

Usual semiclassical arguments lead to the trace formula 

d{E) = d{E) + d°'^{E), (6.1) 

where d{E) is the smooth part of the level density, and d°^'^{E) is its fluctuating part. 

From [3] and [2], it follows that d{E) is given by the following Weyl law (for a circle 
of radius R and A = k in ( \l.2h ): 

f.27r I ^ 

1 



d(E) ^ \ 



d(f) . 



(6.2) 



Vi + /2(0) 

Due to the dependence of A on k, the mean part of level counting function, N{E), is 
not just the integral of d{E) as usual, but is slightly different : 



N{E) ^ -{kRf + —kR 



4tt 



2lT 



d(f) 



l + /(0)^-/(0) 



(6.3) 



The fluctuating part of the level density, d°^'^{E), is given to leading order by the sum 
over all periodic orbits (see e.g. [16]): 



A 



^2^ 



TC 



--{Rp)e 



i[fcip-/ip+7r/4] 



+ C.C. 



(6.4) 



where Ip is the periodic orbit length, Ap is the area swept by the periodic orbit family, fip 
is the phase accumulated due to the caustics, Rp is the total reflection coefficient for a 
given periodic orbit equal the product of reflection coefficients in all points of reflection, 
and (Rp) is its average value over all initial points. 
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Figure 9. Geometrical reflection from the circle. 



Periodic orbits for the circle are regular polygons characterized by two integers 
and M. The integer N gives the number of reflections with the boundary, and the 
integer M determines the number of full rotations around the origin. For co-prime 
and M, the periodic orbit is primitive. Otherwise, it corresponds to the r**^ repetition 
of a primitive periodic orbit where r = (M, A^) is the largest common factor of M and 
N. 

For the circle, 



where 



Ip = 2RN sm9M,N , Ap = tcR^ sin^ 9m,n , /^p = —N (6.5) 

dM,N = . (6.6) 

For the problem under consideration, each reflection with the boundary corresponds to 
the following reflection coefficient [2]: 

R = ^^I^lzl^ (6.7) 
sm + 1 cos (p 

where is the polar angle of the collision point and 6 is the angle between the trajectory 
and the tangent at the point of incidence (see Fig. [9]). 

For a periodic orbit determined by integers M, N, the total reflection coefficient is 
the product of the reflection coefficients for all points of the collisions 

RpW = II -■ / , , o zi V (6-8) 

-L-^ sm 6m,n + 1 cos(0 + 2n6M,N) 

where 9m,n is given by (16.61) . and is the polar angle of the initial incident point. 

For co-prime M and A^, this product is calculated analytically in Appendix D , with 
final result 

f 1 N = 2q 

RpW) = i sinhAAr-i(-l)gcos</.jV = 9^ i 1 " ^^'^^ 

I sinhAAf+i(-l)9cos(/)Ar ^(/ "I" -L 
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Here, A is determined from the equation sinhA = sm6M,N, or 



A = In (^^1 + sin^^M,7v + sin^M,7vJ • (6.10) 

For a repetition of a primitive periodic orbit (i.e. when M = vMq and = tNq with 
(Mo, Nq) = 1), the reflection coefficient is 

Rpr = {Rp{<P)r (6.11) 

where Rp{(f)) is the reflection coefficient for the primitive orbit determined by Mq and 
A^o. 

To calculate the mean reflection coefficient it is necessary to integrate the above 
formulas over all initial incidence angles 0: 

(Rpr) = riRpm^ . (6.12) 



27r 

This integral can be calculated analytically for all r: 

I^RJW'^ = (-1)' + .inh A/V/^.-. (J^) (6.13) 

where P2r-i{x) is a polynomial of degree 2r — 1 of the variable x = 1/ coshAA^. 
In particular, for r = 1,2,3,4, 

Pi{x) = 2x, P3(x) = -4x^ P^ix) = 12x^-8x^+2x, P7(x) = -40x^+48x^-16x1(6.14) 

For large odd A^ and flxed r, the average reflection coefficient is exponentially close to 
the Neumann value 1, 

{Rp) ^ 1 - 4r\-^^^ . (6.15) 

For the triangular periodic orbit with A^ = 3 and M = 1, the average reflection coefficient 
is also very close to 1: {R3) ~ .96391, but for its repetitions it starts to deviate from 
it. For example, for the second, third, and fourth repetition, it has the following values: 
.85972, .69842, and .49720. 

Each term in the oscillating part of the trace formula (16.41) is of the order of 
fc"^/^, which is the dominant contribution when k ^ 00. In general, there exist terms 
decreasing as higher power of k. Usually, these are just small corrections to the existing 
periodic orbit amplitudes and are rarely taken into account. For example, trace formulas 
for odd and even states are slightly different due to mainly the existence of rogue states 
[2]. These eigenfunctions exist only for even states and are described by 

vl/(r, 0) = Mkor) + :M^Ji(A;r) cos0, (6.16) 

where ko is determined from the condition Jiik^R) = 0. Since Jq{x) = — Ji(x), these 
states are the exact analogues of m = states for pure Neumann boundary conditions 
and, as it is easy to check, give corrections proportional to k~^ to the diameter orbit and 
its repetitions corresponding to Ip = 2Rj with integer j. Notice that peaks with odd j 
are absent in the usual trace formula (16. 4p which is typical for desymmetrized systems. 
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In the problem we consider here the situation is different. The shift of the 
whispering gallery mode (15. 4p produces unusual peaks in the trace formula. According 
to (15.41) . the level density of the new levels are given by the formula 

d{x)= Six-xp)= (c+^^^S{P-xC -g{x) + 6), (6.17) 



p 

where 



1 . 1 1 



7r/2 



C=l + -/ = - + -/ Vl + cos2 0d0 . (6.18) 

TT 2 TT Jq 

Using the Poisson summation formula we conclude that 

C - + ^ j (1 + 2 E cos(L(x)r)) (6.19) 

' r=l 

with 

L{x) = 2tx{xC - -^x^/^ ^ 1 ^ _ (g 20) 

Therefore in the limit a; ^ oo these states correspond to a new periodic orbit with the 
length equal to 

/ = 2nRC ^ 6.92i? . (6.21) 

At finite x, the peak associated with this orbit is slightly smaller due to the x^^^ 
correction which dies slowly. 

The numerically-computed length density of orbits is shown in Fig. [TUl For 
comparison, the length density for the Neumann boundary conditions is also represented 
in the same figure but, for clarity, with the opposite sign. All the peaks coincide with 
those of the circle except the additional one associated with the orbit (I6.2ip indicated 
by the arrow. The analysis of the amplitudes of the peaks (e.g. for the triangular orbit 
and its repetitions) confirms fl6.12p . 



7. Conclusion 



In this paper, we studied properties of a circular quantum billiard with specific 
boundary conditions introduced in |2j. Due to the explicit dependence of the boundary 
conditions on momentum, the semiclassical limit of this model is unusual and is neither 
integrable nor chaotic. Following [2], we denote such systems as near integrable. Their 
characteristic property is the strong localization of wave functions in the space of 
azimuthal quantum numbers. The main reason of such behaviour in this system is 
the formal analogy between the recurrence relations for coefficients of the eigenfunction 
expansion, and the one- dimensional Anderson model. In a sense, our system is similar 
to kicked systems where the localization has been established in [13], the role of kicks 
being played by collisions with the boundary. 

The eigenfunction with momentum k can be localized with azimuthal quantum 
number m between and kR. States localized far from the boundaries of this interval 



^ 



7000 



-7000 




Figure 10. Upper curve: oscillating part of the length density. Lower curve: the same 
but for Neumann boundary conditions. The arrow indicates additional orbit (j6.2ip . 



decay (in the mean) exponentially from the point of localization with localization length 
of the order of 1. Eigen- momenta of these states are close to zeros of the derivative of 
the Bessel function J^(^) with m equal the localization center. 

States locahzed close to the boundaries may deviate from pure exponential 
localization. In particular, states having large components with m near kR exponentially 
decrease only for smaller m. The large-m behaviour of these states is described by a 
continuous approximation, and they have oscillations for m in between kR and \/2kR, 
decreasing quickly only for m > \/2kR. The possibility of formation of fractal states 
similar to those investigated in |19] requires additional study. 

Among other consequences of strong localization it is worth mentioning the almost- 
degeneracy of states with different symmetries, and the Poissonian character of spectral 
statistics for energy eigenvalues with the same symmetry. 

As for usual quantum dynamical systems, it is possible to write down a semiclassical 
trace formula relating the quantum spectrum with the sum over periodic orbits. For 
the problem considered, all periodic orbits but one are the same as the integrable case 
of the circular billiard with Neumann boundary conditions, but their amplitudes are 
different due to a different coefficient of reflection with the boundary. The exceptional 
orbit is related to a partially evanescent mode and its length is unusual (see fl6.2ip . 

Though throughout the paper we focus only on a particular example of boundary 
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conditions (with /(0) = cos(0)), our discussion is general and the generahzation for 
other boundary functions should be straightforward. 
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Appendix A. Vanishing of a determinant 

Consider the {2p+ 1) x {2p+ 1) tri-diagonal matrix whose near-diagonal elements equal 
1, and diagonal elements Vn, enumerated consecutively from —p to p, obey 

Vo = 0, (A.l) 

and 

l^_„ = -K for n = l,...,p. (A.2) 

The purpose of this Appendix is to prove that the determinant of this matrix is zero. To 
do this, we show that this matrix always has an eigenvector, a„, with zero eigenvalue. 

By definition a^+i = 0, so it is always possible to find quantities a„ with n = p, . . . ,0 
such that the following relations are fulfilled: 

a„+i + a„_i + \4 a„ = (A. 3) 

Indeed, fixing Op = 1 from these recurrence relations, it follows that ap_i = — V^, 
= Vp-iVp — 1 and so on, until a certain uniquely defined expression for qq (related 

to continued fractions). 

Now consider the same recurrence relation (1A.3P but for negative n = —p, ~{p — 



1), . . . , 0. Because, by the assumption ( ]A.2[) that V^n = —Vn, one can choose solutions 
of (lA.Sp for negative n in the form 

a.n = (-l)"a„ (A.4) 
with above. 

Due to the condition flA.ll) . the remaining equation flA.3|) with n = is also fulfilled. 



proving the existence of an eigenvector with zero eigenvalue and the vanishing of the 
determinant of the matrix considered. 



Appendix B. Effective action integral 

The purpose of this Appendix is to calculate the integral playing the role of the effective 
action for the whispering gallery mode in Section H] for x oo: 

J{x) = I arccos(prrt(a;))dm . (B.l) 

J mi 
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Here, Pmi^) is defined in f l2.3p . mi = mi(x) and m2 = m2{x) are solutions of equations 

pmi(a;) = -l, p„2(x) = 1 . (B.2) 

Tlieir asymptotic values are given by (14.181) and (14.241) correspondingly. 
To find the asymptotics of the integral (]B.1|) . we split it into 3 parts: 

J{x) = Ji{x) + J2{x) + Mx) . (B.3) 

Here Ji{x) denotes the integral over the region close to mi = rup + 1 where Pm{x) can 
well be approximated by the pole term 

Jlix) = I arccos dm (B.4) 

imp+i V rn-rupj 

and Ml — nip ^ 1 . 

J2{x) is the integral over the region far from the pole but close to x where pm{x) is 
well described by the uniform approximation (14. 4p 

Z"*^^ / Ai(y) \ 
J2{x) = / arccos I J dm, (B.5) 



I Ah V vM{y)) 

where y = v = {x/2)^^^, and M2 is chosen such that x^/'^ ^ M2 — x ^x. 

Finally, J3(x) describes the integral over large m where the approximation (14.30 is 
valid: 



m2 I 1^2 3, 



Ma 



J3(x) = / arccos W — - 1 + — — — dm . (B.6) 



x^ 2(m^ — X 



The integral Ji{x) is straightforward, and one finds 

JAx) = (tc — arccos —)du = tt(Mi — nip — 1) (B.7) 

./1 u 



— (Ml — nip) arccos — h In ( Mi — nip + J (Mi — m. 



I — nip * ^ y 

In the limit Mi — m^ ^ 1, we obtain 

Ji(x) = |(Mi - nip) + ln(Mi - mp) + In 2 - vr + 1 . (B.8) 



As f = 


{x/2) 


(IB.5) 






r-M2 


J2{x) = 






J Ah 



dm = -(M2-Mi)+lnAi(^ )-lnAi(^ ) .(B.9) 



V Ai V J 2 V V 
In the second Airy function one can use the approximation 

Ai{y)^Ai'{-7]){y + 7]), (B.IO) 

where as above —rj denotes the first zero of the Airy function, and in the first Airy 
function it is possible to use the asymptotic formula 

Ai(l/)^-=^^y-^/^exp(-^y3/^) . (B.ll) 
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Thus one finds 

J2(x) ^ -(M2-Mi)-ln M -ln[2v/^Ai'(-r/)]-- In ( ) -3 ) -(^ 

In the integral ( 1B.6P the second term is a small correction and 

dm . (B.13) 



arccos \ / — — 1 



x2 2v/2-mVa;2(m2 -x2) 



Substituting m = xt, 

Jz{x) = X 



( rM2/x\ . 1^ 

/ -/ arccosTt^^dt — / , .(B.14) 



After simple transformations we obtain that Jsix) for x 00 has the following 
asymptotics: 

Mx) ^Ix- |(M2 - x) + 1^{M2 - xf'^ + Un {^^^ (B.15) 



where 

/ = ^ arccos Vt^ - 1 &t = V2E (^-^ - | ^ .3393026 (B.16) 
and E{k) is the complete elliptic integral of the third kind, 

For later use we note that 

^ ^ Vl + cos2 0d0 = 72^ (^-^ . (B.18) 

Combining all the terms together, dependence on Mi and M2 indeed disappears as 
it should, and 

J(x) ^ I X + ax^^^ + -\nx + P (B.19) 

6 



where 



a = ^ 2.915016, (B.20) 

5 1 

/3 = -— In 2 - TT - - In TT - In Ai'(-r7) ^ -2.64782 . (B.21) 
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Assume that ?/ is a root of the derivative of the Bessel function Jp{x), 

J'M = . (c.i) 

The purpose of this Appendix is to calculate Jp+k{x) at the point x = y, taking into 
account two terms of the semiclassical expansion when y oo under the assumption 
that p/x is finite but k <^ x. Our starting point is the asymptotic formula for the Bessel 
function, slightly more accurate than (12 .7^ . which can be found e.g. in 



Jp{x) 



(^2_p2)-l/4 



cos ^p{x) 



bi(x) 



\/ x'^ — p' 



sin ^p{x) + 



(C.2) 



where h\ depends only on x/p 
o\[x) 



(C.3) 



8 24(1 -xVp2) ' 

and 0p(x) is the same as in (12.81) . 

Difi^erentiating this expression over one finds that its root, obeys the equation 



sin $p(|/) 



1 



&i(y) 



2(y2-p2)J 



Using the semiclassical expansion of one gets 



COS$p(?/) . 



V 

^ ^p{x) - k arccos - + ■ _ _ . 

X 2^x2 -p2 



(C.4) 



(C.5) 



Then 



Jp+k{x) 



x'-ip+kfr'^' 



cos I $.n(x) — k arccos — h 



X 



+ 



— , sm (PJx) — k arccos — 



2v/a;2-p2 

(C.6) 



Up to the second power of 1/x, this equals 

"2 



Jp+kix) ^ \l -[x~ -p 



2 „2n-1/4 ^ 



pk 



2a/x2 — p2 



COS I $p(x) — k arccos 



+ 



bi 



a/x2 — p2 2A/a;2 — p2 



sin I $p(x) — A; arccos 



• (C.7) 



Expansion of the trigonometric functions and the definition of y ( 1C.4P leads to Eq. (13. 8p . 
Using 2Jp = Jp_i — Jp+i and the recurrence relations for the Chebyshev polynomials 

P,.+i(x) + Pk-i{x) = 2xPkix) (C.8) 

(where Pk{x) stands for the Chebyshev polynomials of the first and the second kind), 
one may also obtain Eq. (13. 9p . 
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To find the total reflection coeflicient for a given periodic orbit characterized by two 
integers M and N, it is necessary to calculate the product (16.81) 

(D.l) 



A J. air 



n=0 



sin 9 



M,N 



+ icos(0 + 26' 



M,N) 



where 6 



M,N 



ttM/N. Denoting z = e^<^, 6 



M,N 



N-1 



-ifR{z) = n 



n=0 



^2g4ien ^ 2iz sin ^e^^^" 



x2e4ien _ 2iz sin 6'e2i^" + 1 ' 



9, and -R(0) = R{z), one flnds 
+ 1 



(D.2) 



Expanding the quadratic polynomials in the numerator and the denominator into 
products of their roots, we obtain 



:-ifR{z) = n 



n=Q 



{Y^-m){Yn + w) 
{Yn + m){Yn-w) 



where Yr, 



ze^'^'\ u 



1/f, and 
a/i + sin^ 9 + sin 6* . 



(D.3) 



(D.4) 



When M and are coprime integers, all products in the above formula can be reduced 
to the following product: 

N-1 



n=Q 



^2-Ki/Nn\ 



with a certain value of t. But the last product equals — 1, therefore 



irR{z) 



{-iz/u) 



N 



[iz/u) 



N 



iiz/v) 



N 



-iz/v) 



N 



Finally 



where 



i?(0) 



A = In 



1, 

sinhAiV-i(-l)(^-i)/^cos0Af 
sinhAAr+i(-l)(^-i)/2cos <I)N ' 



when 
when N 



even 
odd 



1 + sin^ -kM/N + sin ttM/N 



[D.5) 



(D.6) 



(D.7) 
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